A=wk1read('C:\Documents and Settings\jchays\Desktop\Research\Spatial Project\ALM Paper\data\ALM_ML_IRE.wk1',1,0);
W1=wk1read('C:\Documents and Settings\jchays\Desktop\Research\Spatial Project\ALM Paper\data\ALM_W_IRE.wk1');
% Matlab cannot read excel files, but Excel can save data als lotus files with the extension wk1
% 1,1 means skip the first row (e.g. if the file contains variable names
% and the first column (e.g. if the file contains region names)
%
% program to estimate panel data model (N regions*T time periods) 
% data sorted first by time and then by spatial units, so first region 1,
% region 2, et cetera, in the first year, then region 1, region 2, et
% cetera in the second year, and so on
% by ols, regional fixed or random effects, with and without spatial autocorrelation
%
% written by: J.Paul Elhorst 9/2004
% University of Groningen
% Department of Economics
% 9700AV Groningen
% the Netherlands
% j.p.elhorst@eco.rug.nl
%
% Note: some routines or functions should be downloaded from
% James. P LeSage's website www.spatial-econometrics.com
%
% dimensions of the problem
N=15; % number of spatial units
T=12; % number of time periods
K=24; % number of explanatory variables
nobs=N*T;
% normalize W if still necessary
W=normw(W1); % LeSage
y=A(:,[3]); % column number in the data matrix that corresponds to the dependent variable
x=A(:,[4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27]); % column numbers in the data matrix that correspond to the independent variables
idt = eye(T);
fw = kron(idt,W);
% ols estimator
results=ols(y,x); % LeSage
vnames=strvcat('ALM','constant','ALM_1','gdp_pc','ur','density','deind','trade','fdi','oldage','leftc','tcdemc','llvote','govcon','y88','y89','y90','y91','y92','y93','y94','y95','y96','y97','y98'); % should be changed if x is changed
prt_reg(results,vnames); %LeSage
% extension with spatial lag
info.lflag=0;
info.rmin=-1;
info.rmax=0;
results=sar_panel(y,x,W,T,info); % elhorst
prt(results,vnames,1); % LeSage
% ----------------------------------------------------------------------------------------
% fixed effects, within estimator
% demeaning of the y and x variables
meany=zeros(N,1);
meanx=zeros(N,K);
for i=1:N
    ym=zeros(T,1);
    xm=zeros(T,K);
    for t=1:T
        ym(t)=y(i+(t-1)*N,1);
        xm(t,:)=x(i+(t-1)*N,:);
    end
    meany(i)=mean(ym);
    temp1=mean(xm);
    for j=1:K
        meanx(i,j)=temp1(j);
    end
end
ee=ones(T,1);
ywith=y-kron(ee,meany);
% The first explanatory variable is assumed to be the intercept and will be eliminated
xtemp=x(:,[2:K]);
meanxtemp=meanx(:,[2:K]);
xwith=xtemp-kron(ee,meanxtemp);
results=ols(ywith,xwith);
vnames=strvcat('ALM','ALM_1','gdp_pc','ur','density','deind','trade','fdi','oldage','leftc','tcdemc','llvote','govcon','y88','y89','y90','y91','y92','y93','y94','y95','y96','y97','y98'); % should be changed if x is changed
prt_reg(results,vnames);
FE=meany-meanxtemp*results.beta; % including the constant term
yme = y - mean(y);
error=y-kron(ee,FE)-xtemp*results.beta;
rsqr1 = error'*error;
rsqr2 = yme'*yme;
FE_rsqr2 = 1.0 - rsqr1/rsqr2 % r-squared including fixed effects
% ----------------------------------------------------------------------------------------
% fixed effects, within estimator, extension with spatial lag
info.lflag=0;
info.rmin=-1;
info.rmax=0;
results=sar_panel(ywith,xwith,W,T,info);
prt(results,vnames,1);
FE=meany-meanxtemp*results.beta([1:K-1],1)-W*meany*results.rho; % including the constant term
yme = y - mean(y);
for t=1:T
    t1=1+(t-1)*N;t2=t*N;
    error([t1:t2],1)=y([t1:t2],1)-FE-xtemp([t1:t2],:)*results.beta([1:K-1],1)-W*y([t1:t2],1)*results.rho;
end
rsqr1 = error'*error;
rsqr2 = yme'*yme;
FE_rsqr2 = 1.0 - rsqr1/rsqr2 % r-squared including fixed effects
VCVM_TEMP=[results.vcvm];
%VCVM = [VCVM_TEMP(121,121) VCVM_TEMP(12,121) VCVM_TEMP(9,121) ; VCVM_TEMP(121,12) VCVM_TEMP(12,12) VCVM_TEMP(9,12) ; VCVM_TEMP(121,9) VCVM_TEMP(12,9) VCVM_TEMP(9,9)];